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Abstract 



ALHEP is the symbolic algebra program for high-energy physics. It deals with 
amplitudes calculation, matrix element squaring, Wick theorem, dimensional regu- 
f-H \ larization, tensor reduction of loop integrals and simplification of final expressions. 

Qh' The program output includes: Fortran code for differential cross section, Mathe- 

■ matica files to view results and intermediate steps and TeX source for Feynman 

O \ diagrams. The PYTHIA interface is available. 

The project website http://www.riep.by/alriep contains up-to-date executa- 
bles, manual and script examples. 



1 Introduction 
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The analytical calculations in high-energy physics are mostly impossible with- 
out a powerful computing tool. The big variety of packages is commonly used 
[1,2,3,4,5,6,7,8,9,10,11,12,13,14]. Some are general-purpose symbolic algebra 
programs with specific HEP-related plug-ins (REDUCE [15], Mathematica 
[16]), some are designed especially for particle physics (CompHEP [1], SANC 
[2], GRACE [3] etc.) and some are created for specific interaction class or spe- 
cific task. Many of them uses external symbolic algebra core (Form [17], Math- 
Link [16]). They can deal with matrix elements squaring (FeynCalc [7]) or cal- 
culate helicity amplitudes directly (MadGraph [4], CompHEP [1], O'Mega [5]). 
Some packages provide numerical calculations, some require external Monte- 
Carlo generator to be linked. Some programs contain also one-loop calculation 
routines (FormCalc [6], GRACE [3], SANC [2]). Nevertheless, there is no uni- 
form program that meets all the user requirements. 
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Every calculation requires the program-independent check. The optimal tac- 
tics is the simultaneous usage of two (or more) different symbolic algebra 
packages. 

ALHEP is a symbolic algebra program for performing the way from Standard 
Model Lagrangian to amplitude or squared matrix element for the specified 
scattering process. It can also be useful for loop diagrams analysis. The basic 
features are: 

• Diagrams generation using Wick theorem and SM Lagrangian. 

• Amplitude calculation or matrix element squaring. 

• Bondarev functions method for traces calculation. 

• Tensor reduction of loop integrals. 

• Dimensional regularization scheme. 

• Generation of Fortran procedures for numerical analysis (PYTHIA [18] and 
LoopTools [19] interfaces are implemented). 

The current ALHEP version have several implementation restrictions, that will 
be lifted in future. The following features are in progress of implementation: 

• Bremsstrahlung part of radiative correction (the integration over real photon 
phase space). 

• Complete one-loop renormalization scheme including renormalization con- 
stants derivation. 

• Arbitrary Lagrangian assignment. 

After these methods implementation the complete one-loop analysis will be 
available. Please refer to project website for program updates. 

ALHEP website http://www.hep.by/alhep contains the up-to-date executa- 
bles (for both Linux & Win32 platforms), manual and script examples. The 
mirror at http://cern.ch/~makarenko/alhep is also updated. 



2 ALHEP Review. Program Structure 

The ALHEP program internal structure can be outlined as follows: 

• The native symbolic algebra core. 

• Common algebra libraries: 

Dirac matrices, tensor and spinor algebra, field operators & particle wave 
functions zoo. 

• Specific HEP functions and libraries. 

It include Feynman diagrams generation, trace calculations, helicity am- 
plitudes method, HEP-specific simplification procedures, tensor integrals 
reduction and others. 
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• Interfaces to Mathematica, Fortran, TeX and internal I/O format. 

Fortran code is used for further numerical analysis. Mathematica code can 
be used for viewing any symbolic expression in program. But no backward 
interface from Mathematica is currently implemented. TeX output can be 
generated for Feynman diagrams view. Internal I/O format is implemented 
for the most of symbolic expressions, allowing save and restore calculation 
at intermediate steps. 

• Command script processor. 

User interface is implemented in terms of command scripts. The ALHEP 
script language have C-like syntax, variables, arithmetic operations, func- 
tion calls. All HEP-related tasks are implemented as build-in functions. 



2.1 Getting Started 

To use ALHEP one should download the pre-compiled executable for appro- 
priate platform (Linux and Win32 are available) and write a control script 
to describe your task. ALHEP program should be launched with the single 
argument: script file name to be invoked The following steps are required to 
create a workspace: 

• Download ALHEP executables from project website: alhep.gz (for Linux) 
or alhep.zip (Win32). Unpack executable, e.g. 

gzip -d alhep.gz 

• Download up-to date command list: ALHEPCommands.txt. The set of com- 
mands (or options) may be changed in future versions and this manual may 
be somewhat obsolete. Please refer to ALHEPCommands . txt file that always 
corresponds to the latest ALHEP version. The last changes are outlined in 
RecentChanges.txt file at website. 

• Create some working directory and compose command script file therein. 
For example consider the uu — * W + W~ , ~f process with { — h}-helicities of 
initial quarks. To calculate amplitude we create the following script file (call 
it "test.al"): 



SetKinematics (2 , 3 


// 


2->3 process 


,qUARK_U,"p\_l", M e\_l" 


// 


u 


,-QUARK_U, "p\_2" , "e\_2 M 


// 


u-bar 


,WB0Z0N,"f\_l n ,"g\_l" 


// 


w{+} 


, PHOTON, "f\_0 N ,"g\_0" 


// 


photon 


,-WBOZON, "f\_2", "g\_2" ); 


// 


w{-} 


SetDiagramPhysics (PHYS_SM_Q1GEN) 


;// 


SM with 2 quarks only 


SetMassCalcOrder(QUARK_U, 0); 


// 


consider massless 


SetMassCalcOrder(QUARK_D, 0); 


// 


consider massless 


diags = ComposeDiagrams(3) ; 


// 


create diagrams, e~3 


DrawDiagrams (diags , "res.tex", 
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DD.SMALL | DD_SWAP_TALES , FILE.START) ; 

SetFermionHelicity (1 , -1); 

SetFermionHelicity(2, 1); 

SetParameter (PAR_TRACES_BONDEREV, 1) ; 

ampl = CalcAmplitude(RetrieveME(diags) ) ; 

ampl = KinArrange (ampl) ; // arrange result 

ampl = Minimize (ampl) ; // minimize result 

SaveNB ( "res . nb" , ampl , " " , FILE.START I FILE.CLOSE) ; 

f = NewFortranFileC'res.F" , C0DE.F77) ; 

CreateFortranProc(f , "UUWWA", ampl, 
CODE.IS.AMPLITUDE I C0DE_C0MPLEX16 
I CODE_CHECK_DENOMINATORS I CODE_PYTHIA_VECTORS) ; 

The SetKinematics () function declares particles, momenta and polariza- 
tion symbols. The physics is declared with PHYS_SM_Q1GEN option to restrict 
diagrams number. The amplitude for all (d, c, b) internal quarks can be ob- 
tained from generated here by quark mixing matrix replacement: U^ d — > 
^d + ^us + ^ub (chiral limit). Diagrams are created by ComposeDiagramsO 
call. CalcAmplitude () function creates the symbolic value for process am- 
plitude. For detailed discussion of this example see sec. 15.1. 
• Create simple batch file "run. me" like: 
~/alhep_bin_path/alhep test . al 

The ALHEP program creates some console output (current commands, 
scroll bars and some debugging data). If it is not allowed one should redirect 
console output to file here. 

The batch execution with test.al command file takes about 1 minute at 
1.8GHz P4 processor. The following files are created: 

res.nb: Mathematica file containing symbolic expression of amplitude. Cre- 
ated by SaveNB () function. 

res.F: F77 code for numerical analysis created by NewFortranFile () and 
CreateFortranProcO calls. The library file alhep_lib.F is required for 
code compilation and should be downloaded from project website. See sec. 
11.1 for Fortran generation and compilation review. 

res.tex: TeX source for 19 Feynman graphs generated. The AxoDraw [23] 
Latex package is used. The res.tex file should be included into your La- 
TeX document using \input res.tex command. The template document 
to include your diagrams can be found at program website. See fig. 1 in sec. 
15.1 for diagrams generated. 

debug . nb : Mathematica file with debugging information and some interme- 
diate steps. The amount of debugging information is declared in debug . ini 
file in working directory. See sec. 13 for further details. 

See sec. 15 and project website for another examples. It is convenient to use 
some example script as template and modify it for your purposes. 
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2.2 Calculation scheme 



The usual ALHEP script contains several steps: 

(1) Initialization section 

Declaration of process kinematics: initial and final-state particles, titles 
for particle momenta & polarization vectors (SetKinematics ()). 

Physics model definition. The SM physics or part of SM Hamiltonian 
should be specified (SetDiagramPhysicsO). The shorter Hamiltonian is 
selected, the faster is Wick theorem invocation. 

Polarization declaration. Every particle is considered as polarized with 
abstract polarization vector by default. The specific helicity value can be 
set manually or particle can be marked as unpolarized. The several ways 
of polarization involving are inplemented, see sec. 4.4 (SetPolarizedO, 
SetFermionHelicityO, SetPhotonHelicityO, ...). 

Setting mass-order rules for specific particles (SetMassCalcOrder ()). 
One can demand the massless calculation for light particle, that greatly 
saves evaluation time. One can also demand keeping particle mass with 
specific order M n and drop out the higher-order expressions like M n+1 . 
It allows to consider the leading mass contribution without calculating 
precisely. 

(2) Diagrams generation 

Feynman diagrams are generated for specific e n order using the Wick 
theorem algorithm (ComposeDiagramsO) User may draw diagrams here 
(DrawDiagramsO), halt the program (HaltO) and check out if diagrams 
are generated correctly. 

After the diagram set is generated one may cut-off not interesting di- 
agrams to work with the shorter set or select the single diagram to work 
with (SelectDiagramsO). The loop corrections are calculating faster 
when processed by single diagrams. 

Then matrix element is retrieved from diagrams set (RetrieveMEO). 

Before any operation with loop matrix element one should declare 
the N — dimensional space (SetNDimensionSpace ()). Some procedures 
involve N— dimensional mode automatically for loop objects, but most 
functions (arranging, simplification etc.) don't know the nature of expres- 
sion they work with. Therefore the dimensional mode should be forced. 

The diagrams set and any symbolic expression may be saved and re- 
stored at next session to refrain from job repetition (SaveO, LoadO). 

(3) Amplitude calculation 

The amplitude evaluation (CalcAmplitude ()) is the faster way for 
multi-diagram process analysis. 

All the particles are considered polarized. The spinor objects are pro- 
jected to abstract basis spinor. The basis vectors are generated in nu- 
merical code to meet the non-zero denominators condition. See sec. 7 for 
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method details. 

(4) Matrix element squaring (coupling to other) 

The squaring procedure (SquareMEO) is controlled by plenty of op- 
tions, intended mostly for the performance tuning and debugging. 

It basically includes reduction of gamma-matrix sequences coupled to 
kinematically dependent vectors. It reduces the number of matrices in 
every product to minimum. The item-by-item squaring is followed. For 
loop x born* couplings the virtual integrals are involved. See sec. 8 for 
details. 

Amplitude calculation is fast, but its result may be more complicated 
than squaring expressions. Amplitude depends on particle momenta, po- 
larization vectors and additional basis vectors. For unpolarized process 
the averaging cycle is generated in Fortran code, and complex-numbers 
calculation should be performed. The squaring result for unpolarized pro- 
cess is the polynomial of momenta couplings only. Hence there is no 
unique answer what result is simpler for 10-15 diagrams reaction. One 
should definitely use amplitude method for more than 10-15 diagrams 
squaring. 

(5) Loop diagrams analysis 

The tensor virtual integrals are reduced to scalar ones (Evaluate (), 
sec. 9). The scalar coefficients of tensor integral decompositions are also 
reduced to scalar integrals in the most cases (for 1 — 4 point integrals). 

For scalar loop integrals the tabulated values are used. There is no 
reason to tabulate integrals with complicated structure. Hence scalar in- 
tegrals table contain the A and B integrals with different mass configu- 
ration. It also contains a useful Dq chiral decomposition. Other integrals 
should to be resolved using LoopTools-like [19,20] numerical programs. 

The renormalization procedure is under construction now. The counter- 
terms (CT) part of Lagrangian leads to CT diagrams generating. In the 
nearest future the abstract renormalization constants (S m , 5f etc.) will be 
involved and tabulated for the minimal on-shell scheme. The automatic 
derivation of constants is supposed to be implemented further. Please 
refer ALHEP website for implementation progress. 

(6) Simplification 

The kinematic simplification procedure is available (KinSimplifyO, 
sec. 12). It reduces expression using all the possible kinematic relations 
between momenta and invariants. The minimization of +/* operations 
in huge expressions can also be performed (Minimize ()). 

(7) Fortran procedure creation for numerical analysis 

F90 or F77- syntaxes for generated procedures are used. Generated code 
can be linked to PYTHIA, LoopTools and any Monte-Carlo generator for 
numerical analysis. 



6 



3 ALHEP script language 



The script syntax is similar to C/Java languages. Command line breaks at 
";" symbol only. Comments are marked as // or /*. . .*/. Operands may be 
variables or function calls. If no function with some title is defined, the operand 
is considered as variable. The notation is case-sensitive. 

The script language have no user-defined functions, classes or loop operators. 
It seems to be useless in current version. 

All ALHEP features are implemented as build-in functions. The execution 
starts from the first line and finishes at the end of file (or at Halt () command). 

Variable types are casted and checked automatically (in run-time mode). No 
manual type specifying or casting are available. List of ALHEP script internal 
types: 

• Abstract Symbolic Expression (expr). 

Result of function operations. Can be stored to file(SaveO) and loaded 
back(LoadO). Basic operations: +,-,*. The division is implemented using 
Frac(a,b) function. Are not supposed to be inputed manually, although a 
few commands for manual input exist. 

• Integer Value (int). 

Any number parameter will be casted to integer value. Basic operations: 
+ ,-,*,|. The division is performed using Frac(a,b) function. No fractional 
values are currently available. 

• String (str). 

String parameters are started and closed by double quotes ("string vari- 
able"). Used to specify symbol notation (e.g. momenta titles), file names, 
etc. Basic operations: -(-(concatenate strings). 

• Set of Feynman Diagrams (diagrams) . 

Result of diagrams composing function. Basic operations: SaveO, LoadO 
and SelectDiagramsO. Diagrams can also be converted into TeX graphics. 

• Matrix Element (me). 

Contains symbolic expression for matrix element and information on it's 
use: list of virtual momenta (integration list) etc. 

The total up-to-date functions list can be found in ALHEPCommands . txt file 
at project website. 
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4 Initialization section 



4-1 Particles 

Particles are determined by integer number, called particle kind (PK). The 
following integer constants are defined: ELECTRON, MUON, TAULEPTON, PHOTON, 
ZBOZON, WBOZON, QUARK.D, QUARK.U, QUARK.S, QUARK.C, QUARK.B, QUARK.T, 
NEUTRINO.ELECTRON, NEUTRINO.MUON, NEUTRINO.TAU. The ghost and scalar 
particles are not supposed to be external and their codes are unavailable. 
Antiparticles have negative PK that is obtained by " -PK" operation. 

If kinematics is declared the particles can be secelted using the particle ID 
(PID) number. Initial particles have negeative ID (-1, -2) and final are pointed 
by positive numbers (1, 2, 3...). 



4-2 Kinematic selection 

Before any computations may be performed one needs to declare the kinematic 
conditions. They are: number of initial and final state particles, PK codes, 
symbols for momentum and polarization vector for every particle. 



SetKinematics( (int)N_Initials , (int)N_Finals , 

[(int)PK_I, (str)momentum_I , (str)polarization_I , ...]); 

Here N_Initials and N_Finals - numbers of initial and final particles in kine- 
matics. The next parameters are particle kind, momentum and polarization 
symbols, repeated for every particle. For example, the e~(A; 1 , ei)e + (k 2 , e 2 ) — > 
fj>~(jPi, e 3 )/i + (j»2, ^a) process should be declared as follows: 

SetKinematics (2 , 2, 

ELECTRON, "k\_l" , "e\_l" , 
-ELECTRON, "k\_2" , "e\_2" , 
MUON, "p\_l", "e\_3", 
-MUON, "p\_2", "e\_4"); 



4-3 Particle masses 

All the particles are considered massive by default. Particle can be declared 
massless using the SetMassCalcOrder function. 
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Set MassCalcOrder ( (int ) PK , (int ) order) ; 
PK: particle kind, 

order : maximum order of particle mass to be kept in calculations. Zero value 
means massless calculations for specified particle. Negative value declares 
the mass-exact operations. 

The mass symbols are generated automatically and look like m e , mw etc. 
Hence, all the electrons in process will have the same mass symbol. To declare 
unique mass symbols for specific particles the SetMassSymO function is used. 
The different masses for unique particles are often required to involve the 
Breit-Wigner distribution for particles masses. 



SetMassSym((int)PID, (str or expr)mass); 

PID: particle ID in kinematics (< for initial and > for final particles), 
mass: new mass symbol, like "m\_X" for m X - 

4-4 Polarization data 

All the particles are considered as polarized initially. The default polarization 
vector symbols are set together with kinematic data (in SetKinematicsO). 

To declare particle unpolarized the SetPolarizedQ function is used: 



SetPolarized((int)PID, (int)polarized); 



PID: particle ID in kinematics, 

polarized: 1 (polarized) or (unpolarized). 

One can declare the specific polarization state (helicity value) for particles. 

The photon (k,e) helicities are involved in terms of two outer physical mo- 
menta: 

e t = ((P+ • k )P-» - (P- ■ k)p+^ ± ie lia p v p%jP.k v ') /N Kp± . (1) 



SetPhotonHelicity((int)PID, (int)h, (str)"pP", (str)"pM"]) ; 
PID: photon particle ID in kinematics, 

h: helicity value: ±1 or 0. Zero value clears the helicity information 
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pP, pM: p± base vectors in (1) formula. The 7-coupled forms of (1) (precise 
or chiral [26]) are used automatically if available. 

One may select the transverse unpolarized photon density matrix (instead of 
usual -(l/2)g tlv ): 



P h b P b b b 2 



SetPhotonDMBase((int)PID, (str)" P" ) ; 



PID: photon particle ID in kinematics, 
"P" : Basis vector P in (2). 



The fermion (k,e) helicity is declared using: 



SetFermionHelicity((int)PID, 


(int)h); 


SetFermionHelicity((int)PID, 


(srt)"h"); 



PID: fermion particle ID in kinematics, 

h: helicity value: ±1 or 0. Zero value clears the helicity information. Density 

matrix is usual: ^Zuu — > fry± (massless fermion). 
"h" : symbol for scalar parameter in the following density matrix (massless 

fermion): J2 u u (l/2)p(l + ^75) • 

Notes: 

• SetPolarized(PID, 0) call will also clear helicity data. 

• SetXXXHelicity(PID, 1, ...) also sets particle as polarized (previous 
SetPolarized(PID,0) call is canceled). 

• SetXXXHelicity(PID, 0, . . . ) clears helicity data but does not set parti- 
cle unpolarized. 



4-5 Physics selection 



One can either declare the full Standard Model (in unitary or Feynman gauge) 
or select the parts of SM Lagrangian to be used. The QED physics is used by 
default. The Feynman rules corresponds to [27] paper. 



SetDiagramPhysics( (int)physID ); 
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physID: physics descriptor, to be constructed from the following flags: 
PHYS.QED: Pure QED interactions, 
PHYS.Z and PHYS.W: Z- and W-boson vertices, 

PHYS.MU, PHYS.TAU: Muons and tau leptons (and corresponding neutrinos), 
PHYS.SCALARS: Scalar particles including Higgs bosons, 
PHYS.GHOSTS: Faddeev-Popov ghosts, 
PHYS_GAUGE_UNITARY: Use unitary gauge, 

PHYS.QUARKS, PHYS _ QUARKS _ 2 GEN and PHYS_QUARKS_3GEN: {d, u}, {d, u, s, c} 

and {d, u, s, c, b, t} sets of quarks, 
PHYS_RARE_ VERI CES : vertices with 3 or more SCALAR/GHOST tales, 
PHYS_CT: Renormalization counter-terms (implementation in progress), 
PHYS.SM: full SM physics in unitary gauge (all the flags above except for 

PHYS.CT), 

PHYS_SM_Q1GEN: the SM physics in unitary gauge with only first generation 

of quarks (d, u), 
PHYS_ELW: SM in Feynman gauge with no rare 3-scalar vertice, 
PHYS.EONLY: no muons, tau- leptons and adjacent neutrinos, 
PHYS.NOQUARKS: no quarks, 

PHYS_4B0S0NS_AN0MAL0US: anomalous quartic gauge boson interactions (see 
[28]), affects on AAWW and AZWW vertices. 

The less items are selected in Lagrangian, the faster diagram generation is 
performed. 



5 Bondarev functions 

The Bondarev method of trace calculation is implemented according to [25] 
paper. The trace of 7^-matrices product becomes much shorter in terms of 
F-functions. The number of items for Tr[(l — 7 5 )aia2 • • -a^n] occurs 2 n . For 
example, the 12 matrices trace contains 10395 items usually while the new 
method leads to 64-items sum. 

The 8 complex functions are introduced: 





= 2[(aq^)(bq + )- 


(ae+)(6e_)], 


F 2 


[a,b) 


= 2[(ae + )(6g_)- 


(ag_)(6e+)], 




[a,b) 


= 2[(aq + )(bq_)- 


(ae_)(&e + )], 


F, 


[a,b) 


= 2[(ae_)(bq + )- 


(ag+)(6e_)], 


F 5 


[a,b) 


= 2[(aq_)(bq + )- 


(ae_)(6e+)], 


F 6 


[a,b) 


= 2[(ae_)(6g_)- 


(ag_)(6e_)], 


F 7 


[a,b) 


= 2[(aq + )(bq_)- 


(ae+)(6e_)], 


F s 


(a, 6) 


= 2[(ae + )(bq + )- 


(aq+)(be+)] 



(3) 
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The basis vectors q± and e± are selected as follows: 



ql = i=(l, ±1, 0, 0), e£ = -J=(0, 0, 1, ±i). (4) 



The results for traces evaluation looks as follows: 
Tr[(l - 7 5 )aia 2 ] = Fi(ai, a 2 ) + F 3 (ai, a 2 ), 

Tr[(l - 7 5 )a 1 a 2 a 3 a 4 ] = F^ai, a 2 )Fi(a 3 , a 4 ) + F 2 (a 1; a 2 )F 4 (a 3 , a 4 ) + 
+F 3 (ai, a 2 )F 3 (a 3 , a 4 ) + F 4 (ai, a 2 )F 2 (a 3 , a 4 ). 



Please refer to [25] paper for method details. 



SetParameter (PAR_TRACES_BONDEREV, (int)par) ; 



par : 1 or - allow or forbid Bondarev functions usage. 

The numerical code for F-functions (3) is contained in alhep_lib.F library 
file. The code is available for scalar couplings only: F llv p^q v . If vector Bondarev 
functions remain in result, the Fortran-generation procedure fails. One should 
repeat the whole calculation without Bondarev functions in that case. 



6 Diagrams generation 



The diagrams are generated after the kinematics and physics model are de- 
clared. The Wick-theorem-based method is implemented. 

The only distinct from the usual Feynman diagrams is the following: the vertex 
rules have additional 1/2 factors for every identical lines pair. It makes two 
effects: 

- The crossing diagrams are usually involved if they have different topology 
from original. I.e. if two external photon lines starts from single vertex, they 
are not crossed. Nevertheless all the similar external lines have crossings in 
ALHEP. 

- Some ALHEP diagrams have 2™ factors due to identical internal lines. For 
example, the W + W~ — > { r yyWW vertex} — > 77 — > {'yyWW vertex} — > 
W + W~ diagram will have additional factor 4. The result remains correct 
due to 1/2 factors at every r )^W + W~' vertex. 

The diagrams are generated without any crossings. The crossed graphs are 
added automatically during the squaring or amplitude calculation procedures. 
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diagrams ComposeDiagrams((int)n); 
n: order of diagrams to be created, M x ~*y ~ e n . 

Uses current physics and kinematics information. Returns generated diagrams 
set. 

One may select specified diagrams into another diagrams set: 

diagrams SelectDiagrams( (diagrams)d, (int)iO [, il, i2...]); 

d: initial diagrams set. Remains unaffected during the procedure. 
iO. . iN: numbers of diagrams to be selected. First diagram is "0". 

The new diagrams set is returned. 

To retrieve matrix element from the diagram set the RetrieveMEO is used, 
me RetrieveME( (diagrams)d ); 

d: diagrams list. 

7 Helicity amplitudes 

The amplitudes are calculated according to [24] paper. Every spinor in matrix 
element is projected to common abstract spinor: 



The value of vector Q is selected arbitrary in Fortran numerical procedure. 
The additional basis vector Eq (if exists) is selected meet the polarization 
requirements ((Eq.Q) = 0, (Eq.Eq) = -1). The fractions like l/Tr[P Q P p ] 
may turn 1/0 at some Q and Eq values. The denominators check procedures 
are generated in Fortran code (the C0DE_CHECK_DEN0MINAT0RS key should be 
used in CreateFortranProcO call). If the | Tr[PQP p ] |> 5 check is failed, the 
another Q and Eq values are generated. 




UQU p U p _ iC U Q P p _ U p U p UQ _ iC 

— O . r ___ w i \ i /<-» i ^*T) — — ^ 



(Tr[P Q P p ])V*- 



The projection operator Pq = u p u p is choosen as follows: 



• Pq = Q(l + 75)/2 - for massive external fermions, 

• Pq = Q(l + Eq)/2 - if one of fermions is massless. 
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expr CalcAmplitude( (me)ME ); 
ME: matrix element retrieved from the whole diagrams set. 

The result expression is a function of all the particle helicity vectors. The 
averaging over polarization vectors is performed numerically. The numerical 
averaging procedure is automatically generated in Fortran output if unpolar- 
izaed particles are declared. 

The CODE_IS_AMPLITUDE key in CreateFortranProcO procedure declares 
expression as amplitude and leads to proper numerical code (Ampl x Amlp*). 



8 Matrix Element squaring 

The squaring procedure have the following steps: 

• matrix elements simplification to minimize the 7-matrices number, 

• denominators caching to make procedure faster, 

• item-by-item squaring (coupling to other conjugated), 

• saving memory mechanism to avoid huge sums arranging (SQR_SAVE_MEMORY 
option), 

• virtual integrals reconstruction. 



expr SquareME((me)MEl, [(me)ME2,] [(int)flags ]); 
ME1 : matrix element #1, 

ME2: matrix element #2 (should be omitted for squaring), 
flags: method options (defauls is "0"): 
SQR_CMS: c.m.s. consideration (initial momenta are collinear). The additional 
pseudo-covariant relations (pi.£2 - ► 0,p2-£i - ► 0) appear that simplify work 
with abstract polarization vectors. 
SQR_N0_CR0SSING_1: do not involve crossings for ME #1, 
SQR_N0_CR0SSING_2: do not involve crossings for ME #2, 
SQR.MANDELSTAMS: allow Mandelstam variable usage, 

SQR_PH_GAMMA_CHIRAL: tries to involve photons helicities in short chiral e ± 

form (see sec. 4.4, [26]). 
SQR_PH_GAMMA_PRECISE: involve precise e ± form for photon helicities. If the 

polarization vectors are not coupled to 7^, the vector form is used (see 

sec. 4.4). 

SQR_SAVE_MEMORY: save memory and processor time for huge matrix ele- 
ment squaring. Minimizes sub-results by every 1000 items and skips final 
arranging of the whole sum. No huge sums occurs in calculation in this 
mode, but the result is also not minimal. If result expression contains a 
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sum of 10 5 — 10 6 items (when expanded), the arranging time is significant 
and SQR_SAVE_MEMORY flag should be involved. If no results are calculated 
in reasonable time the CalcAmplitude () (see sec. 7) procedure should be 
used. 

If two matrix elements are given, the first will be conjugated, i.e. the result is 
ME1* x ME2. 



9 Virtual integrals operations 

The tensor virtual integrals are reduced to scalar ones using two methods. 

If tensor integral is coupled to external momentum and p-vector can 

be decomposed by integral vector parameters, the fast reduction is involved. 
The D x integrals for 2 — > 2 process contain the whole basis of 4-dimension 
space and /^-couplings to any external momentum can be decomposed. It 
works well if all the polarization vectors are constructed in terms of external 
momenta. 

The common tensor reduction scheme is involved elsewhere. Tensor integrals 
are decomposed by the vector basis like 



The linear system for scalar coefficients is composed and solved. Implemented 
for Bj and C^j integrals only. The other scalar coefficients should be calculated 
numerically [19]. 



expr Evaluate( (expr)src ); 



Evaluate: Reduction of tensor virtual integrals to scalar ones. The scalar 
coefficients in tensor VI decomposition are also evaluated (not for all integrals). 



expr ConvertInvariantVI( (expr)src ); 



ConvertlnvariantVI: Vector parameters are substituted by scalars. I.e.: 



Note: A- and B- integrals are converted automatically during arrangement. 



expr CalcScalarVI( (expr)src ); 




C (ki, k 2 , mi, m 2 , m 3 ) -> C (kl, (h - k 2 ) 2 , k\, m 1 , m 2 , m 3 ). 
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CalcScalarVI: Substitutes known scalar integrals with its values. The most 
of UV-divergent integrals (A ,B Q ) are substituted. The chiral decomposition 
for _D -integral is also applied. The complicated integrals should be calculated 
numerically [19]. 

The source expressions are unaffected in all the functions above. 



10 Regularization 

The dimensional regularization scheme is implemented. One may change the 
space-time dimension before every operation: 



SetNDimensionSpace( 


(int)val ); 






val: 1 (n-dimensions) 


or (4-dimensions space) 


expr Singular Arrange ( 


(expr)src ); 







Singular Arrange: turns expression to 4-dimensional form. Calculates (n — 4) 1 
factor in every item and drops out all the neglecting contributions. 

SetDRMassPK( (int)PK ) 
SetDRMassPK: set the particle to be used as DR mass regulator. 
PK: particle kind (see sec. 4.1). "0" value declares the default "/i" DR mass. 

11 ALHEP interfaces 

11.1 Fortran numerical code 

The numerical analysis in particle physics is commonly performed using the 
Fortran programming language. Hence, we should provide the Fortran code to 
meet the variety of existing Monte-Carlo generators. 

To start a new Fortran file the NewFortranFileO function is used: 
int NewFortranFile( (str)fn [, (int)type]); 
fn: output file name 
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type : Fortran compiler conventions: C0DE_F77 or C0DE_F90. The FORTRAN 77 

conventions are presumed by default. 
Function returns the ID of created file. 

Then we can add a function to FORTAN file: 



CreateFortranProc( (int)fTD, (str)name, (expr)src [, (int)keys]); 



fID: file ID returned by NewFortranFile () call, 
name : Fortran function name. 

src : symbolic expression to be calculated. The source may be | M T, da/dY 
(use CODE_IS_DIFF_CS flag) or amplitude (CODE_IS_AMPLITUDE flag). The 
result is always da/dY calculation procedure, 
keys: option flags for code generation (default is C0DE_REAL8): 

C0DE_REAL8: mean all the symbols in expression as REAL*8 values. 

C0DE_C0MPLEX16: declare variables type as COMPLEX* 16. 

CODE_IS_DIFF_CS: the source expression is differential cross section. 

CODE_IS_AMPLITUDE the source expression is amplitude and a squaring code 
should be generated: AMPL*DC0N JG (AMPL) . 

CODE_CHECK_DENOMINATORS check denominators for zero. Used to re-generate 
free basis vectors in amplitude code (see sec. 7). 

CODE.LOOPTOOLS: Create LoopTools [19] calls for virtual integrals. 

CODE_SEPARATE_VI create unique title for every virtual integral function (ac- 
cording to parameter values). Do not use with C0DE_L00PT00LS. 

CODE_PYTHIA_VECTORS retrieve vector values from PYTHIA [18] user process 
PUP (I, J) array. 

CODE_POWER_PARAMS factorize and pre-calculate powers if possible. 
C0DE_N0_4VEC_ENTRY do not create a 4-vector entry for function. 
C0DE_N0_C0NSTANTS do not use predefined physics constants. All variables 
becomes external parameters. 
C0DE_N0_SPLIT do not split functions by 100 lines. 
C0DE_N0_C0MM0N: don't use CONNON-block to keep internal variables. 

The scalar vector couplings, Bondarev functions (see sec. 5) and SabcdVilfyp'i'Pt 
objects are replaced with scalar parameters to be calculated once. These func- 
tions are calculated using alhep_lib.F library procedures. 

Some compilers works extremely slow with long procedures. Therefore Fortran 
functions are automatically splitted after every 100 lines for faster compilation. 

The internal functions functions have TEMPXXXO notation. To avoid problems 
with several ALHEP-genarated files linking we should rename the TEMP-prefix 
to make internal functions unique. 



SetParameterS(PAR_STR_FORTRAN_TEMP, (str)prefix); 
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prefix: "TMP1" or another unique prefix for temporary functions name. 



To obtain the better performance of numerical calculations ALHEP provides 
the mechanism for minimization of +/x operations in the expression. We 
recommend to invoke the Minimize () function before Fortran generation (see 
sec. 12). 

11.1.1 PYTHIA interface 

The PYTHIA [18] interface is implemented in terms of UPINIT/UPEVNT proce- 
dures and 2 — > iV phase-space generator. 

The momenta of external particles are retrieved from PYTHIA user process 
event common block (PUP (I, J)). The order of particles in kinematics should 
meet the order of generated particles, and the CODE_PYTHIA_VECTORS option 
should be used. 

The template UPINIT/UPEVNT procedures for ALHEP — > PYTHIA junction are 
found at ALHEP website. However one should modify them by adjusting the 
generated particles sequence, including user cutting rules, using symmetries 
for calculation similar processes in single function etc. 

The plane 2 — > N phase-space generator in alhep_lib . F library file is written 
by V. Mossolov. It is desirable to replace the plain phase-space generator by 
the adaptive one for multiparticle production process. 

The more automation will be implemented in future. Please refer to ALHEP 
website for details. 

11.2 Mathematica 

The output interface to Mathematica Notebook [16] file is basically used to 
view the expressions in the convenient form. Implemented for all the symbolic 
objects in ALHEP. 

SaveNB( (str)fn, (expr or me)val [, (str)comm ][, (int)flags ]); 
MarkNB( (str)fn [, (str)comm ][, (int)flags ]); 

fn: Mathematica output file name, 
val : expression to be stored. 

comm: comment text to appear in output file (" " if no comments are required), 
flags: file open flags (default is "0"): 
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0: append to existing file and do not close it afterward, 

FILE_START: delete previous, start new file and add Mathematica header, 

FILE_CL0SE: add closing Mathematica block. 

The MarkNB function is used to add comments only. 

The Mathematica program can open valid files only. The valid x . nb file should 
be started (FILE_START) and closed once (FILE_CL0SE). It is convenient to in- 
sert MarkNB ("x.nb" , " " ,FILE_START) and MarkNB ("x.nb" , " " ,FILE_CL0SE) 
calls to start and end of your script file. 

No backward interface (Mathematica — > ALHEP) is currently available. 



11.3 LaTeX 



The LaTeX interface in ALHEP is only implemented for Feynman diagrams 
drawing. The diagrams are illustrated in terms of AxoDraw [23] package. 

DrawDiagrams( (diagrams)d, (str)fn [, (int)flags] [, (int)draw]); 
MarkTeX( (str)fn, [, (str)comm][, (int)flags]); 

d: diagrams set. 

fn: output LaTeX file name. 

flags: file open flags (O(default): append, FILE_START: truncate old): 
comm: comment text to appear in output file (" " if no comments are required), 
draw: draw options (default is DD.SMALL): 

DD.SMALL: diagrams with small font captions. 

DD_LARGE: diagrams with large font captions. 

DD_M0MENTA: print particles momenta. 

DD_SWAP_TALES: allow arbitrary order for final-state lines (option is included 

automatically for 2 — > N kinematics). 
DD_DONT_SWAP_TALES: deny the DD_SWAP_TALES option. The order of final 

lines is same to kinematics declaration. 



11.4 ALHEP native save/load operations 

The native ALHEP serialization format is XML-structured. It may be edited 
outside the ALHEP for some debug purposes. 
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Save( (str)fn, (diagrams or expr)val); 
object Load( (str)fn ); 

fn: output XML file name. 

val : diagrams set or symbolic expression to be stored. 



12 Common algebra utilities 

The following set of common symbolic operations is available: 

• ExpandO: expands all the brackets and arranges result. Works slowly with 
huge expressions. 

• Arrange (): arranges expression (makes alphabetic order in commutative se- 
quences). The most of ALHEP functions performs arranging automatically 
and there is no need to call Arrange () directly. 

• Minimize (): reduces the number of "sum-multiply" operations in expres- 
sion. Should be used for simplification of big sums before numerical calcu- 
lations. 

• Factor (): factorize expression. 

• KinArrange () : arranges expression using kinematic relations, 

• KinSimplif y () : simplify expression using the kinematic relations. Works 
very slowly with large expressions. Mostly useful for 2 — > 2 process. 



expr Expand ( (expr)src ); 

expr Arrange( (expr)src ); 

expr Minimize( (expr)src [, (int)flags]); 

expr Factor( (expr)src [, (int)flags]); 

expr KinArrange( (expr)src [, (int)flags]); 

expr KinSimplify( (expr) src [, (int)flags]); 

src: source expression (remains unaffected). 
Minimize () function flags (default is MIN.DEFAULT): 

MIN.DEFAULT = MIN.FUNCTIONS I MIN.DENS I MIN.NUMERATORS, 

MIN.FUNCTIONS: factorize functions, 

MIN_DENS: factorize denominators, 

MIN.NUMERATORS: factorize numerators, 

MIN.NUMBERS: factorize numbers, 

MIN_ALL_DENOMINATORS: factorize all denominators I /(a + b) + l/x 

(x + a + b)/(x*{a + b))), 
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MIN_ALL_SINGLE_DENS: factorize single denominators (but not sums, prod- 
ucts etc.), 

MIN_VERIFY: verify result (self-check: expand result back and compare to 
source). 

Factor () function flags (default is 0): 

FACT_NO_NUMBERS: do not factorize numbers, 
FACT_N0_DENS: do not factorize fraction denominators, 
FACT_ALL_DENS: factorize all denominators, 

FACT_ALL_SINGLE_DENS: factorize all single denominators (but not sums, 
products etc.), 

FACT.VERIFY: verify result (self-check: expand result back and compare to 
source). 

KinArrangeO function flags (default is 0): 

KA_MASS_EXACT: do not truncate masses (neglecting SetMassCalcOrder() 
settings), 

KA.MANDELSTAMS: involve Mandelstam variables (for 2^2 kinematics), 
KA_N0_EXPAND: do not expand source. 
KinSimplifyO function flags (default is KS_FACTORIZE_DEFAULT) : 

KS_FACTORIZE_DEFAULT: factorize functions and denominators (the first 
two flags below), 

KS_FACT0RIZE_FUNCTI0NS: factorize functions before simplification, 

KS_FACTORIZE_DENS: factorize denominators (including partial factoriza- 
tion) before simplification, 

KS_FACTORIZE_ALL_DENS: factorize all denominators, 

KS_MASS_EXACT: do not truncate masses (neglecting SetMassCalcOrder() 
settings), 

KS.MANDELSTAMS: involve Mandelstam variables (for 2^2 kinematics), 
KS_N0_EXPAND: do not expand source (if no simplification are found). 



13 Debugging tools 

ALHEP allows user to control the most of internal calculation flow. The debug 
info is stored to debug. nb file (critical messages will also appear in console). 
The debug . ini file contains numerical criteria for messages to be logged, the 
debug levels for different internal classes. Warning: raising the debug . ini val- 
ues leads to sufficient performance drop and enormous debug . nb file growth. 

If one feels some problems with ALHEP usage, please contact author for assis- 
tance (attaching your script file). Do not waste the time for manual debugging 
using debug, ini. 

There are some specific commands to view internal data. For example, the 
whole list of tensor virtual integrals reduction results is kept in internal data 
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storage and can be dumped to Mathematica file for viewing: 

• str ViewParticleData( (int)PID) returns the brief information on parti- 
cle settings. Puts the information string to console and returns it also as a 
string variable. 

PID: particle ID in kinematics. 

• ViewFeynmanRules ( (str)nb_f ile, (int)flags) stores Feynman rules of 
current physics to Mathematica file. 

• ViewTensorVITable( (str)nb_f ile, (int)f lags) ) stores tensor integrals 
reduction table to Mathematica file. The VI reduction table is filled during 
the Evaluate () operation. 

• ViewScalarVTCache( (str)nb_f ile, (int)flags)) stores scalar loop in- 
tegrals values cache. The scalar VI cache is filled during the CalcScalarVI () 
invocation. 

nb_file: output Mathematica file name, 

flags : access flags to Mathematica file: 0,FILE_START and/or FILE.CLOSE. 
Call without parameters turns output to debug . nb file. 



14 System commands 

The two system commands are useful: 

• Halt(): stop further script processing. May be used to test the first part of 
script and save (SaveO) internal result. It first part finished successfully, it 
may be commented (\* . . . *\) and followed by loading procedure (LoadO). 
Then script execution is restarted. 

• Timer (): view time elapsed since the last Timer() call (from program start 
for first call). 



15 Examples 

15.1 Amplitude for qq — > W + W~^ 

Let's consider the example from sec. 2.1 in details. We also extend it for anoma- 
lous quartic gauge boson interactions [28]. And we don't use the Bondarev 
method for traces calculation this time. Please refer to sec. 2.1 for ALHEP 
installation notes. 

We start test . al script from output files creation: 
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nbfile = "uuWWA_MPXXX.nb" ; // Mathenatica file name 

MarkNB(nbfile, "", FILE.START) ; // Create file 

texfile = "res.tex"; // LaTeX file name 

MarkTeX (texfile, "", FILE.START); // Create file 

Then 2 — * 3 process kinematics and physics are declared: 



u(pi, ei) «(p2, e 2 ) -> 7(/ , g ) W + (f u W (/ 2 , # 2 ). (5) 

SetKinematics (2 , 3 // 2->3 process 

,QUARK_U,"p\_l","e\_l M // u 

,-QUARK_U,"p\_2", M e\_2" // u-bar 

,WB0Z0N,"f\_l","g\_i" // W{+} 

,PH0T0N,"f\_0","g\_0" // photon 

,-WBOZON, "f\_2", "g\_2" ); // W{-} 
SetDiagramPhysics (PHYS_SM_Q1GEN I PHYS_4B0S0NS_AN0MAL0US) ; 

We declare physics with u- and <i-quarks only. The amplitude will be sum- 
marized for all the possible internal quarks numerically. It requires the simple 
replacing of quark mixing matrix in resulting Fortran code: U^ d — > U^ d + U^ s + 

u 2 ub ). 

Next we declare the u- and rf-quarks massless: 

SetMassCalcOrder (QUARK_U, 0); // consider massless 

SetMassCalcOrder (QUARK_D, 0); // consider massless 

Set polarizations to " -+UUU" : 

SetFermionHelicity (1 , -1); // u 

SetFermionHelicity(2, 1); // u-bar 

Create diagrams set and store it to LaTeX file: 

diags = ComposeDiagrams(3) ; //e~3 order 

DrawDiagrams (diags , texfile); 

Next we include the following lines: 

Save ( "diags .xml" , diags) ; //save to XML file 

//HaltO; //stop execution 

//diags = Load( "diags .xml") ; //load from XML file 

We can save diagrams, stop the program now and view diagrams generated. 
To stop ALHEP session the HaltO line should be uncommented. Then we 
modify our script as follows: 
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Fig. 1. The diagrams generated for uu — > W W +r y process (see res.tex file). The 
anomalous quartic gauge boson interaction affects the first two diagrams. 

/* diags = ComposeDiagrams(3) ; // commented 

... // commented 

Halt(); */ // commented 

diags = Load (" diags. xml" ) ; // uncommented 

If we run the script again, it will skip the diagrams generation step and load 
diagrams from XML file. 

Matrix element retrieval: 

me = RetrieveME (diags) ; //get matrix element 

SaveNB(nbf ile, me, "Matrix element"); //view 
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Calculate helicity amplitude, arrange result and minimize the +/x operations 
number: 



ampl = CalcAmplitude (me) ; 

SaveNB(nbf ile, ampl, "Amplitude after CalcAmplitude ()") ; 
ampl = KinArrange(ampl) ; 

SaveNB(nbf ile, ampl, "Amplitude after KinArrangeO ") ; 
ampl = Minimize (ampl) ; 

SaveNB(nbf ile, ampl, "Amplitude after MinimizeO ") ; 

The another breakpoint can be inserted here. The result for amplitude is saved, 
the Halt and Load commands are commented for further use: 

Save (" ampl. xml " , ampl); // save amplitude 

//MarkNB(nbfile, FILE.CLOSE) ; Halt(); // close NB and exit 
//ampl = Load( "ampl. xml" ) ; // load amplitude 

This breakpoint allows to repeat the next Fortran creation step without re- 
calculating of matrix element. 

Let's average over final state polarizations in further numerical procedure. Set 
final particles unpolarized: 



SetPolarized(-l, 0) 
SetPolarized(-2, 0) 
SetPolarized(-3, 0) 



// set unpolarized 
// set unpolarized 
// set unpolarized 



The Fortran output for differential cross section: 

SetParameterS (PAR_STR_FORTRAN_TEMP , "TMP1 11 ) ; 

f = NewFortranFileC'uuWWA.F", C0DE_F77) ; 1 1 ±11 file 

CreateFortranProc(f , "uuWWA" , ampl, 

CODE_IS_AMPLITUDE| //square amplitude 

C0DE_CHECK_DEN0MINAT0RS I //check 1/0 limits 
C0DE_C0MPLEX16 I //complex values 

C0DE_P0WER_PARAMS I //F(M~2) instead of F(M) 

C0DE_PYTHIA_VECT0RS) ; //use PYTHIA PUP (I, J) vectors 

The SetParameterS call sets the unique notation for internal variables and 
functions. Please do not make it too long. The complex-type code is required 
for proper amplitude calculation. 

Close Mathematica output file at the end of script: 

MarkNB(nbfile, FILE.CLOSE) ; 

The execution of this script takes less than 2 minutes at 1.8GHz P4 processor. 
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We will not discuss the structure of generated uuWWA.F file in details. But 
some remarks should be done: 

Line 5: The main function call. The following parameters are declared (order 
is changed here): All the parameters (except the kQOrig) are of COMPLEX* 16 
type. Ones the C0DE_C0MPLEX16 option is set, all the real objects are treated 
as complex. 

kQOrig (INTEGER) : The ID of w(first)-quark in PYTHIA PUP(I,J) array. 

Possible values: 1 or 2. 
PAR_a_0, PAR_a_c, PAR_a_n, PAR_ah_c, PAR_ahat_n: Anomalous quartic 

gauge boson interaction constants ao, a c , a n , a c , a n [28]. 
PAR_CapitalLambda: Scale factor A for anomalous interaction [28]. 
PAR_VudP2 : Quark mixing matrix element squared |t/ MC ;| 2 . The U% d + U£ s + 
U^ b value may be passed to summarize the whole diagrams (neglecting 
quarks masses). The numbers for mixing matrix elements may be obtained 
using QMIX_VAL(ID1,ID2), QMIX_SQR_SUM(ID) and QMIX_PROD_SUM func- 
tions of alhep_lib.F library. 
Line 29: Internal COMMON-block with PAR (XX) array. All the scalar cou- 
plings and other compound objects are precalculated and stored in PAR (XX) . 
Lines 49-55: External momenta initialization from PYTHIA PUP (I, J) ar- 
ray. The order of external vectors is expected as follows: 
PUP (1,1) ,PUP(I,2) : initial particles. If kQ0rig=2 the order is backward: 

PUP (I, 2), PUP(I.l). 
PUP (1, 3.. N): final particles in the same order as in SetKinematics () 
call. One should modify this section (or SetKinematics () parameters) 
to make the proper particles order. 
DO 10 1=1,4 

p_l(I) = DCMPLX(PUP(I,kQ10rig) ,0D0) // kQlOrig = kQOrig 
f_0(I) = DCMPLX(PUP(I,4) ,0D0) 

p_2(I) = DCMPLX(PUP(I,kQ20rig) ,0D0) // kQ20rig = 3-kQ0rig 
f_l(I) = DCMPLX(PUP(I,3) ,0D0) 
f_2(I) = DCMPLX(PUP(I,5) ,0D0) 
10 CONTINUE 

Lines 69-231: Polarization averaging and basis vector generation cycle. For 
any momenta set the PAR (XXX) array is filled. Then denominator checks 
and amplitude averaging are performed. 

Lines 244-252, 438-444, ... Interaction constants and particle masses defi- 
nitions in sub-procedures. The constants can be declared as main function 
parameters using C0DE_N0_C0NSTANTS option in CreateFortranProc func- 
tion. 

For the complete pp — > W + W^ / y analysis the following steps are required: 

• The PYTHIA client program should be written. The template files are avail- 
able at ALHEP website. 
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The another helicity configuration +-UUU should be calculated separately. 
The another channels q^j — > W + W~^ (i ^ j) should be calculated and 
included into generator. 



15.2 Z- boxes for e e + — > ji ji + 

Let's calculate some box diagrams now. Consider the following process: 

e~(pi,ei) e + (p 2 ,e 2 ) -> fT(fi,gi) // + (/2,# 2 )- (6) 

As in previous example, we start command script from files initialization: 

nbfile = "Zbox.nb"; // Mathenatica file name 

MarkNBCnbfile, "", FILE.START); // Create file 

texfile = "res.tex"; // LaTeX file name 

MarkTeX(texfile, "", FILE.START); // Create file 

The 2^2 kinematics declaration: 

SetKinematics(2, 2, // 2->2 process 

ELECTRON, "p\_l", "e\_l" , // e~{-} 

-ELECTRON, "p\_2 M , "e\_2", // e~{+} 

MUON, "f\_l", "g\_l", // imr{-} 

-MUON, "f\_2", "g\_2"); // mu~{+} 

The Standard model physics in Feynman gauge (we omit quarks for faster 
diagrams generation): 

SetDiagramPhysics (PHYS.ELW I PHYS.NOQUARKS) ; 

Set leptons massless and declare the iV-dimensional space: 

SetMassCalcOrder (ELECTRON, 0); //massless electrons 

SetMassCalcOrder (MUON, 0); // massless muons 

SetNDimensionSpace (1) ; 

Use Mandelstam variables throughout the calculation: 
SetParameter (PAR.MANDELSTAMS, 1) ; 
Consider unpolarized particles: 

SetPolarizedd , 0); // unpolarized e~{-} 

SetPolarized(2, 0); // unpolarized e~{+} 
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SetPolarized(-l, 0) ; 
SetPolarized(-2, 0) ; 



// unpolarized mu"{-} 
// unpolarized mu"{+} 



Compose born and one-loop diagrams: 

diags_born = ComposeDiagrams (2) ; //e~2 order 

DrawDiagrams (diags_born, texf ile) ; 

diags_loop = ComposeDiagrams (4) ; //e~4 order 

Save("diags_loop.xml" ,diags_loop) ; 

/ / diags_loop=Load ( " diags_loop . xml " ) ; 

DrawDiagrams (diags_loop, texf ile) ; 

The 220 loop diagrams are created, saved to internal format and TeX file. 
The ComposeDiagrams (4) procedure takes 5-10 minutes here. It is convenient 
to comment ComposeDiagrams (4) -Save () lines at second run and Load loop 
diagrams from disk. 




Fig. 2. Born level diagrams for e e + — ► fi /i 




Fig. 3. Part of 220 loop diagrams stored to res .tex file 

Let's select the double Z-exchange box graphs from the whole set (194 and 
195 diagrams at fig. 3): 



diag_box = SelectDiagrams (diags_loop , 194, 195) ; 

Next we couple the loop and born matrix elements: 

me_born = RetrieveME(diags_born) ; 
me_box = RetrieveME(diag_box) ; 
me_sqr = SquareME(me_box, me_born) ; 
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The simplification procedures are not included into SquareME implementation. 
It may take much time to arrange items in huge expression. Therefore all the 
simplification procedures are optional and should be called manually: 

me_sqr = KinArrange(me_sqr) ; 
me_sqr = KinSimplif y (me_sqr) ; 

SaveNB(nbf ile, me_sqr, "squared & simplified"); 

The reduction of tensor virtual integrals follows: 

me_sqr = Evaluate (me_sqr) ; 
me_sqr = KinArrange(me_sqr) ; 
me_sqr = KinSimplif y(me_sqr) ; 
SaveNB(nbf ile, me_sqr, "VI evaluated"); 

Next we convert scalar integrals to invariant-dependent form and replace with 
tabulated values: 

me_sqr = ConvertlnvariantVI (me_sqr) ; 

me_sqr = CalcScalarVI (me_sqr) ; // use pre-calculated values 
me_sqr = KinArrange(me_sqr) ; 
SaveNB(nbf ile, me_sqr, "VI scalars "); 

Turn to 4-dimensional space, drop out (n — 4) 1 items and final simplification: 

me_sqr = SingularArrange (me_sqr) ; 
SetNDimensionSpace (0) ; 
me_sqr = KinArrange(me_sqr) ; 
me_sqr = KinSimplif y(me_sqr) ; 

Save result and create Fortran code with LoopTools [19] interface: 

SaveC'ZBox.xml" ,me_sqr) ; // save result 

//me_sqr = LoadC'ZBox.xml") ; // reload result 

SaveNB(nbf ile, me_sqr, "Z boxes result"); // view result 
f = NewFortranFileC'ZBOX.F", C0DE_F77) ; 
CreateFortranProc(f , "ZBOX", me_sqr, 

C0DE_P0WER_PARAMS I C0DE_L00PT00LS) ; 

View tensor integrals reduction table and close Mathematica output file: 

ViewTensorVITable(nbf ile) ; 
MarkNBCnbfile, FILE.CLOSE) ; 

The script runs about 15 minutes on 1.8GHz P4 processor. The half of this 
time takes the ComposeDiagrams(4) procedure. 
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Fig. 4. The result expression in Zbox.nb file. This result contains UV-regulator 
term x, that should cancel one in Bq(s, Mz, Mz) integral. It can be checked using 
GetUVTermO function. 

Code remarks: 

Line 5: The main function call. The 3 parameters are usual Mandelstam 
variables (s,t,u, type is complex). The current ALHEP version does not 
care about interdependent parameters in Fortran output. And all the three 
Mandelstam variables may occur in parameters list. The future versions will 
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be saved from this trouble. 
Line 14,50: Include LoopTools header file ("looptools.h"). See LoopTools 

manual [19] for details. 
Line 20,21: Retrieve LoopTools values for UV-regulator getdeltaO and 

DR-mass squared getmudimO . 

The complete code of examples including scripts, batches and output files are 
available at ALHEP website. 



16 Conclusions 



The new program for symbolic computations in high-energy physics is pre- 
sented. In spite of several restrictions remained in current version, it can be 
useful for computation of observables in particle collision experiments. It con- 
cerns both multiparticle production amplitudes and loop diagrams analysis. 

The nearest projects are: 

• Bondarev functions method improvement, 

• Complete renormalization scheme for SM, 

• Complete covariant analysis of the one-loop radiative corrections including 
the hard bremsstrahlung scattering contribution. 

• Arbitrary Lagrangian assignment. 

Refer ALHEP project websites for program updates: 

http://www.hep.by/alhep , 

http : / /cern . ch/~makarenko/ alhep . 
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